function [a0, a2, a4] = fn_fit_poly(vx, vep)

n_right = find(vep>0, 1);
n_left = n_right - 1;

vind = [(n_left - 4):n_left, n_right:(n_right + 4)];

x = vx(vind);
y = vep(vind);

Num = length(x);

V0 = ones(1,Num);
V2 = x.^2;
V4 = x.^4;

n = @(Vec1, Vec2) Vec1*transpose(Vec2);

M = [n(V0,V0), n(V0,V2), n(V0,V4); ...
    n(V2,V0), n(V2,V2), n(V2,V4); ...
    n(V4,V0), n(V4,V2), n(V4,V4)];
V = [n(V0,y); n(V2,y); n(V4,y)];

tv = M\V;
a0 = tv(1);
a2 = tv(2);
a4 = tv(3);

end